1 Find EICs that look like contaminants

library(xcms)
library(MetabolomiQCsR)
library(purrr)
library(dplyr)
library(magrittr)
library(chemhelper)
library(massageR)
library(RColorBrewer)
library(heatmaply)

F first we read a raw file

file <- "Z:/_Data/LIP1/0036_SevaMeal/mzML/0036_LIP1p_20160201_034_2_1_1.mzML"
raw <- xcmsRaw(file, profstep = 0)

raw
## An "xcmsRaw" object with 2530 mass spectra
## 
## Time range: 3.1-1260 seconds (0.1-21 minutes)
## Mass range: 99.9961-1698.7933 m/z
## Intensity range: 200-12561400 
## 
## MSn data on  0  mass(es)
##  with  0  MSn spectra
## Profile method: bin 
## Profile step: no profile data
## 
## Memory usage: 169 MB

Then we can look for contaminant-like features

contaminants <- EIC_contaminants(raw, min_int = 1E4)
## Allowing multi-threading with up to 8 threads.
contaminants
## # A tibble: 341 × 2
##          mz                  EIC
##       <dbl>               <list>
## 1  100.0753 <tibble [2,530 × 3]>
## 2  102.1272 <tibble [2,530 × 3]>
## 3  109.1009 <tibble [2,530 × 3]>
## 4  111.0909 <tibble [2,530 × 3]>
## 5  113.0690 <tibble [2,530 × 3]>
## 6  113.1070 <tibble [2,530 × 3]>
## 7  118.0861 <tibble [2,530 × 3]>
## 8  118.1219 <tibble [2,530 × 3]>
## 9  119.0852 <tibble [2,530 × 3]>
## 10 121.0509 <tibble [2,530 × 3]>
## # ... with 331 more rows

2 Plot contaminant EICs

plot_chrom(contaminants$EIC[[20]], RT_col = "scan_rt", Intensity_col = "intensity")

3 Annotate contaminants with known contaminants

We can now match the list of proposed contaminants to a list of known common contaminants.

cont_list <- get_cont_list("positive") %>% extract2(1)

cont_names <- db.comp.assign(   mz = contaminants$mz,
                                rt = rep(0,nrow(contaminants)),
                                comp_name_db = with(cont_list, paste0(`Compound ID or species`," (",`Ion type`,")")  ),
                                mz_db = cont_list$`Monoisotopic ion mass (singly charged)`,
                                rt_db = rep(0,nrow(cont_list)),
                                mzabs=0.01,ppm=15,
                                ret_tol=Inf
                            ) %>% 
              data_frame(name = .)
##                                         [,1]
## Unique hits                               22
## Non-unique hits                            4
## Compounds not found                      735
## Markers had multiple compounds assigned    0
contaminants %<>% bind_cols(cont_names)
 
contaminants
## # A tibble: 341 × 3
##          mz                  EIC         name
##       <dbl>               <list>        <chr>
## 1  100.0753 <tibble [2,530 × 3]> NMP ([M+H]+)
## 2  102.1272 <tibble [2,530 × 3]> TEA ([M+H]+)
## 3  109.1009 <tibble [2,530 × 3]>             
## 4  111.0909 <tibble [2,530 × 3]>             
## 5  113.0690 <tibble [2,530 × 3]>             
## 6  113.1070 <tibble [2,530 × 3]>             
## 7  118.0861 <tibble [2,530 × 3]>             
## 8  118.1219 <tibble [2,530 × 3]>             
## 9  119.0852 <tibble [2,530 × 3]>             
## 10 121.0509 <tibble [2,530 × 3]>             
## # ... with 331 more rows

4 Heatmap of contaminants

4.1 Calculate EIC correlation matrix

corr <- map(contaminants$EIC, ~ .x$intensity) %>% do.call(cbind,.) %>% cor

4.2 Calculate dendrograms

z <- heat.clust(corr,
                scaledim =    "none", 
                zlim     =    c(-Inf,Inf),
                zlim_select = "",
                reorder  =    c("column","row"),
                scalefun =    function(x) x,
                distfun =     function(x) dist(x, method="euclidean"),
                hclustfun =   function(x) hclust(x, method="complete")
)

4.3 Make heatmap

# row/colnames need to be unique
colnames(z$data) <- contaminants %>%    mutate(name = paste0(name,", ",round(mz,4))) %>% 
                                        extract2("name") %>% 
                                        make.unique(sep="_")

rownames(z$data) <- colnames(z$data)


heatmaply(z$data,
                   Rowv=z$Rowv,
                   Colv=z$Colv,
                   symm=TRUE,
                   scale="none",
                   colors=rev(colorRampPalette(brewer.pal(10, "RdBu"))(40)),
                   column_text_angle = 45
          ) %>% 
layout(margin = list(l = 300, b = 150)) # fix cut of labels